FERMILAB-PUB-93/066-A 
March 1993 



Axion Miniclusters and Bose Stars 

Edward W. Kolb«>( 2 ) and Igor I. Tkachev«<( 3 ) 



W NASA/Fermilab Astrophysics Center 
Fermi National Accelerator Laboratory, Batavia, IL 60510 

^Department of Astronomy and Astrophysics, Enrico Fermi Institute 
The University of Chicago, Chicago, IL 60637 

^Institute for Nuclear Research of the Academy of Sciences of Russia, Moscow 

117312, Russia 



Evolution of inhomogeneities in the axion field around the QCD 
epoch is studied numerically, including for the first time important non- 
linear effects. It is found that perturbations on scales corresponding to 
causally disconnected regions at T ~ 1 GeV can lead to very dense axion 
clumps, with present density p a > 10~ 8 gcm~ 3 . This is high enough for 
the collisional 2a — > 2a process to lead to Bose-Einstein relaxation in the 
gravitationally bound clumps of axions, forming Bose stars. 
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The invisible axion is one of the best motivated candidates for cosmic dark mat- 
ter, despite being subject to strong cosmological and astrophysical constraints on its 
properties (10 10 GeV < f a < 10 12 GeV for the axion decay constant; 10~ 5 eV < m a < 
10~ 3 eV for the axion mass) [T[]. As dark matter, axions would play a role in the 
evolution of primordial density fluctuations and formation of large scale structure. In 
addition to its generic properties, axions also have unique features as dark matter. 
For instance, large amplitude density fluctuations produced on scales of the horizon 
at the QCD epoch [@] lead to tiny gravitationally bound "miniclusters" ||. It was 
found that the density in miniclusters exceeds by ten orders of magnitude the local 
dark matter density in the Solar neighborhood 0. This might have a number of 
astrophysical consequences, as well as implications for laboratory axion searches [|J. 

In previous studies, spatial gradients of the axion field in the equations of motion 
were neglected. This is a reasonable assumption for temperatures below the QCD 
scale where the evolution of coherent axion oscillations can be treated as pressureless, 
cold dust. However, we find that just at the crucial time when the inverse mass of the 
axion is approximately the size of the horizon, gradient terms become important, and 
a full field-theoretical approach is needed. Here we present the results of a numerical 
study of the evolution of the inhomogeneous axion field around the QCD epoch. 
Though we only consider spherically symmetric configurations, the importance of 
the combined effect of the field gradients and the non-liner attractive self interaction 
should also occur if we relax spherical symmetry. The resulting axion clumps are 
much denser than previously thought, reaching the critical conditions for Bose star 
formation ||. 

The axion field 8(x) is created during the Peccei-Quinn symmetry breaking phase 
transition at T ~ /„, uncorrelated on scales larger than the horizon at this time ||. 
For T < f a , the field becomes smooth on scales up to the horizon, H~ l (T), where 
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H is the expansion rate. This continues until T = T\ « 1 GeV when the axion mass 
switches on, i.e., when m a (T\) « 3H(Ti). Coherent axion oscillations then transform 
fluctuations in the initial amplitude into fluctuations in the axion density. 

Since the initial amplitude of coherent axion oscillations on the horizon scale 
-f/" _1 (Xi) is uncorrelated, one expects typical positive density fluctuations on this 
scale will satisfy p a ~ 2p a , where p a is mean cosmo logical density of axions ||. At 
the temperature of equal matter and radiation energy density, T e = 5.5 Q a h 2 eV [0, 
these fluctuations are already non-linear and will separate out as miniclusters with 
p a ~ 3(10eV) 4 ~ 10~ 14 gcm -3 0. The minicluster mass will be of the order of 
the dark-matter mass within the Hubble length at temperature T\, M mc ~ 1O _9 M . 
The radius of the cluster is R mc ~ 10 13 cm, and the gravitational binding energy will 
result in an escape velocity v e /c ~ 10~ 8 . Note that the mean phase-space density 
of axions in such a gravitational well is enormous: n ~ p a m~ 4 v~ 3 ~ 10 48 / 4 2 ; where 
fi2 = f a /l0 12 GeV. 

We will show below that due to non-linear effects, a substantial number of regions 
can have axion density at T > T e many times larger than 2 p a . 

Let us parametrize the energy density of a single fluctuation as p a (T e < T < 
Ti, 6j) = 3 $(6 l j)T e s/4, where 6^ is the misalignment angle at T 1; s is the entropy den- 
sity, and = 1 corresponds to the mean axion density. The energy density inside 
a given fluctuation is equal to the radiation energy density at T = $(0j)T e . At that 
time the fluctuation becomes gravitationally non-linear and collapses. Consequently, 
at T e 

Pa (d t ) ~ $ 4 (^)p a (T e ), (l) 

will be the minicluster density after it separates out as a bound object. Even a 
relatively small increase in is important because the density depends upon the 
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fourth power of 

Ref. || demonstrated that due to anharmonic effects for fluctuations with 6i close 
to 7T, some correlated regions can have values of larger than just a factor of two. 
The reason is simple: the closer 8i is to the top of the axion potential, 

V{6) = m 2 a (T)f a (l - cosfl) = A*(T)(1 - cosfl), (2) 

the later axion oscillations commence. However this effect alone is not very significant. 
In the range 0.1 < £ < 10 -3 we can parametrize it as $(6^) 1.5(#j/7r) 2 £~ - 35 , where 
£ = (tt — 9i)/ir, and is significantly larger than 2 only for field values very finely 
tuned to the top of the potential. Moreover, the axion field is not exactly coherent 
on the horizon scale, and small fluctuations might spoil this picture. 

At temperatures T ^> T%, the potential is negligible in the equations of motion 
compared to the gradient terms which force the field to be homogeneous on scales 
less than the horizon. At T -C T\, on the contrary, gradients can be neglected and 
one can treat the evolution of fluctuations as that of a pressureless gas. Clearly at 
T ~ Ti, both the gradient terms and the potential are important, and in order to 
find the energy density profile at freeze out one has to trace the inhomogeneous field 
evolution through the epoch T ~ Tj. 

It is convenient to work in conformal coordinates with metric ds 2 = a 2 (r])(dr) 2 — 
dx 2 ). During radiation dominance a oc 77 and 77 oc T~ l . The dependence of the axion 
mass upon the temperature at T > Aqcd can be found in the dilute-instanton-gas 
approximation ||, and can be parametrized as a power law, m 2 (r]) = m 2 (ri*)(ri/i]*) n , 
where n = 7.4 ± 0.2 0. Introducing the field ip — l^i the equations of motion for 
a spherically symmetric axion fluctuation in an expanding Universe is of the form 
ip — ip" — 2ip'/ r + f] n+3 sin (ip/rj) = 0, where fj is the reduced conformal time parameter 
V = v/v*! an d i^aijl*) = H{ri*). The radial coordinate r is defined in the comoving 
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FIG. 1. Energy density contrast in a fluctuation with initial radius r = 1.8 and 
Of = 2.75 at several moments of time as a function of comoving radius r. The 
density contrast is normalized to the value of the homogeneous energy density at 
fj = 4. 

reference frame, with r = 1 corresponding to R p hy S {v*) = H ^iv*)- 

We integrated this equation numerically for a wide range of initial conditions. 
We evolved configurations which are at rest at fj = 0.1. The initial distribution of 
the field can be parametrized by the initial radius of the fluctuation, ro, the initial 
value of the field inside, Of, the initial value of the field outside, Of", and the width 
of transient region, Ar. The important common feature is that the final density 
distribution develops a sharp peak in the center. The larger the gradients of initial 
configuration, the higher the final peak, e.g., the peak grows with increase in \0f — 0f\. 



4 



FIG. 2. Energy density profiles at fj = 4 for identical initial fluctuations evolved with 
different Lagrangians. Solid line: axion case; dashed line: V(9) oc 9 2 /2; dotted line: 
V{9) oc 9 2 /2 + 9 A /4; dash-dotted line: axion potential with field gradients switched 
off. 

The peak also grows with decreasing width of the transient region. We present here 
the results of runs with initial amplitude of the field outside the fluctuation equal to 
the r.m.s. value of the misalignment angle, i.e., 9f = tt/V3, and width of transient 
layer Ar ~ 0.6. 

Energy density profiles as a function of time are presented in Fig. 1 for a typical 
case. At fj — 1 there are two waves, incoming and outgoing, both propagating with 
the velocity of light. At approximately fj — 2 the incoming wave reaches the center 
and the outgoing wave reaches r 3.5. At later times the wave front does not move 
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significantly because the axion mass effectively switches on at fj m 2, and the edge of 
the fluctuation "freezes." 

One reason for energy density growth at later times is the continuing increase of 
the axion mass. However the relative density contrast in the center with respect to 
the unperturbed homogeneous environment continues to increase up to the final time 
of integration, fj = 4. This is entirely a non-linear effect. One can see this in the 
following way: The average pressure over a period of homogeneous axion oscillations 
in potential Eq. (||) is negative, and is equal to (P) ~ — A^(T)#q/64, where 9 is 
the amplitude of the oscillations 0. In other words, the axion self- interaction is 
attractive. The larger the amplitude of oscillations inside the fluctuation, the more 
negative will be the pressure inside, and consequently, fluctuations with excess axions 
will contract in the comoving volume. In addition, matter with a smaller pressure 
suffers less redshift in the energy density. To see this effect we present in Fig. 2 the 
final density profiles correspondin g to identical initial field distributions evolved with 
different potentials: the axion potential of Eq. (Q); the axion potential with gradients 
artificially switched off; a pure harmonic potential, V(6) oc 6 2 /2, where (P) = 0; 
and the potential V(8) oc 8 2 /2 + 8 4 /A, where (P) > 0. Note that for the harmonic 
potential, at fj = 4 the maximal density excess is only about 3, i.e., ten times smaller 
than for the axion potential. 

The dependence of the energy density contrast in the center upon the initial radius 
is shown on Fig. 3. In the whole range of values of r$ plotted, the energy density takes 
its maximum value just in the center of the final configuration. Only if r$ < 1.55 
or r > 2.05 does the final energy density profile have a maximum at some non- 
zero radius. In a sense, the initial radius of the fluctuation in the plotted range is 
more or less tuned in such a way that the arrival of the incoming wave at the center 
is synchronized with the switching on of the axion mass. However, there is nothing 
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FIG. 3. Dependence of density contrast in the center of a fluctuation at fj — 4 upon 
the initial radius of a fluctuation for several values of the initial misalignment angle 
inside the fluctuation. 

unnatural in this "synchronization," since as larger and larger scales enter the horizon 
in an expanding Universe there will always be a scale for which the incident wave of 
a disappearing fluctuation reaches the center just at the moment of freeze out. 

Quantitatively, the assumption of spherical symmetry is very important. How- 
ever, in general any isolated contrast in the initial misalignment angle will decay via 
incoming and outgoing waves which will not possess spherical symmetry. The overall 
picture will be the same as in the spherical case, but the values of the maximal energy 
contrast in the final configuration at a given Of will be smaller. Note in this respect, 
that the final density contrast rapidly grows with increase of Of (see Fig. 3) due to the 
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attractive self-interaction resulting in negative pressure. This has nothing to do with 
the symmetry of the fluctuation, and we may expect to find large density contrasts 
in regions where the field values happen to be close to tt initially [ |10| |. 

The effect of the field gradients is important not only in the discussion of the 
formation of high density peaks, but also in the careful estimate of the mean density 
of axion matter. We found that the total excess mass of axions within a fluctuation, 
compared to the homogeneous background, does not vary much, and is equal approx- 
imately to half of the excess mass if gradients in the equations of motion would be 
neglected. This deficit might be attributed to the redshift at early times, fj < 2, when 
axions are still relativistic. 

The energy density contrast plotted in Figs. 1 and 2 will coincide with the factor 
§(6i) in Eq. (||) if we assume that the mean cosmological density of axions corresponds 
to homogeneous oscillations with initial amplitude equal to the r.m.s. value of the 
misalignment angle. As we have noted already, the energy density in an axion clump 
after it separates out from the general expansion will be $ 4 (#j) times larger than the 
energy density at T e . So a density contrast of 30 will correspond to roughly a factor 
of 10 6 in the energy density of the cluster at T < T e . 

All axion miniclusters could be, in principle, relevant to laboratory axion search 
experiments, since for a minicluster with $ as small as 2, the density is 10 10 larger than 
the local galactic halo density. However, the probability of a direct encounter with 
a clump is small. The interesting question arises, could there be any astrophysical 
consequences of very dense axion clumps? Below we shall discuss the possibility of 
"Bose star" formation inside axion miniclusters. 

The physical radius of an axion clump at T e is larger by many orders of magnitude 
than the de Broglie wavelength of an axion in the corresponding gravitational well. 
Consequently, gravitational collapse of the axion clump and subsequent virialization 
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can be described in the usual terms of cold dark matter particles. In a few crossing 
times some equilibrium (presumably close to an isothermal) distribution of axions 
in phase space will be established. It is remarkable that in spite of the apparent 
smallness of axion quartic self-couplings, |A | = {fn/ fa) 4 ~ lO -53 /-^ 4 , the subsequent 
relaxation in an axion minicluster due to 2a — > 2a scattering can be significant as a 
consequence of the huge mean phase-space density of axions ||. In the case of Bose- 
Einstein statistics the inverse relaxation time is (1 + n) times the classical expression, 
or r^ 1 ~ nv e ap a /m a , where a is the corresponding cross section. For part icles 
bounded in a gravitational well, it is convenient to rewrite this expression in the form 



The shallower the gravitational well at a given density of axions, the larger the mean 
phase space density, and consequently the smaller the relaxation time due to the v\ 
dependence in Eq. (|3|) . Note also the dependence of the inverse relaxation time upon 
the square of the particle density. 

The relaxation time (|3|) is smaller then the present age of the Universe if the 
energy density in the minicluster satisfies 



where pio = p/(10eV) 4 and V-$ = v e /10~ 8 . If this occurs, then an even denser core 
in the center of the axion cloud should start to form. An analogous process is the 
so-called gravithermal instability caused by gravitational scattering. This was studied 
in detail for star clusters, where the "particles" obey classical Maxwell-Boltzmann 
statistics. Axions will obey Bose-Einstein statistics, with equilibrium phase-space 
density n(p) = n con( j + [e l3E — l] -1 , containing a sum of two contributions, a Bose con- 
densate and a thermal distribution. The maximal energy density that non-condensed 
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7\-2 -2 2 
A Pa V e . 



(3) 




(4) 
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axions can saturate is p t her ~ m a v li which corresponds to n t her ~ 1- Consequently, 
given the initial condition n ^> 1, one expects that eventually the number of particles 
in the condensate will be comparable to the total number of particles in the region 
where relaxation is efficient. Under the influence of self-g ravity, a Bose star PJTT| 
then forms J5J. One can consider a Bose star as coherent axion field in a gravitational 
well, generally with non-zero angular momentum in an excited energy state ||. 

Comparing Eqs. ([!]) and we conclude that the relaxation time is smaller than 
the present age of the Universe and conditions for Bose star formation can be reached 
in miniclusters with density contrast > 30 at the QCD epoch. For examples of 

such density contrasts, see Figs. 1 and 3. 

Under appropriate conditions stimulated decays of axions to two photons in a 
dense axion Bose star are possible f9,12] (see also 0]), which can lead to the formation 
of unique radio sources — axionic masers. In view of results of present paper we 
conclude that the questions of axion Bose star formation, structure and possible 
astrophysical signatures deserves detailed study. 

In conclusion, we have presented a numerical study of the evolution of inhomo- 
geneties in the axion field around the QCD epoch, including for the first time impor- 
tant non-linear effects. We found that the non-linear effects can lead to a much larger 
core density of axions in miniclusters than previously estimated. The increase in the 
density may be sufficiently large that axion miniclusters might exceed the critical 
density necessary for them to relax to form Bose stars. 

It is a pleasure to thank H. Feldman, J. Frieman, A. Kashlinsky, A. Klypin, D. 
Pogosyan, A. Stebbins, M. Turner, and R. Watkins for useful discussions. This work 
was supported in part by the DOE and NASA grant NAGW-2381 at Fermilab. 
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